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A generalized mean field theory of coarse-graining. 
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Abstract 

A general mean field theory is presented for the construction of equilibrium coarse grained 
models. Inverse methods that reconstruct microscopic models from low resolution experimental 
data can be derived as particular implementations of this theory. The theory also applies to the 
opposite problem of reduction, where relevant information is extracted from available equilibrium 



1-^ ' ensemble data. These problems are central to the construction of coarse grained representations 



of complex systems, and commonly used coarse graining methods are derived as particular cases 
of the general theory. 
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Coarse-graining(CG) methods are widely used in both simulations and theory to con- 
struct simplified models that render complex problems analytically and computationally 
;ractable. Successful examples are the theory of polymers [l|, the theory of liquid crystals 
2| and hydrodynamics |3|. In all these models the system is not described with an atomistic 
resolution, but only a few relevant properties of the system are retained to account for the 
essential physics required to describe the system. 

Two general strategies are commonly used in the construction of CG models. The first in- 
volves finding a solution to an inverse problem, where some macroscopic quantity is known 
and a microscopic structure of the system is required that is consistent with this quan- 
tity. An example of this strategy is the reverse montecarlo (RMC) [J, |5| method where 
a radial distribution function (RDF) is provided as the target property that needs to be 
reproduced. A second class of methods involve reducing the number of degrees of freedom 
present in the system. An example of such an approach is the multiscale coarse-graining 
(MS-CG) method|6|, where a CG model is built from an atomistic trajectory generated from 
a molecular dynamics (MD) simulation. Inverse problems are of primary importance for the 
interpretation of experimental data, where macroscopic quantities are measured and a micro- 
scopic description of the system is required, y, |7| On the other hand, an important challenge 
in data analysis is to develop systematic methods by which the reduction and organization 
of high resolution data can be performed to extract the relevant pieces of information from 
available data. J8|, l9| Extracting relevant quantities is also essential for constructing simpli- 
fied theoretical models that retain the essential aspects of a complex physical process while 
discarding non-relevant degrees of freedom. This is particularly true for complex systems. 



such as biomolecules, where CG models help in exploring large scale properties. [10| 

Several techniques are available in the published literature that deal with the above de- 
scribed problems, but a general formulation that can be used to construct such methods 
does not currently exist. This work establishes a general statistical mechanical formulation 
that can be used to build such models, and which also serves as a general statistical me- 
chanical framework within which the physics of such modeling strategies can be interpreted 
and understood. The theory presented in this work offers a simple and general framework 
for a diverse class of methods which deal with the type of problems described above and is 
intimately connected with traditional classical mean field theories. Thus, the results from 
this work enable the application of techniques from statistical mechanics to develop novel 



and powerful methods for coarse-graining and the study of inverse problems in equilibrium 
statistical mechanics. 

A primary goal in the construction of a CG model is to provide a microscopic description 
of the system under examination which can reproduce the ensemble average of a given set 
of macroscopic observables, ({o"j}). An average of such a variable in the Gibbs-Boltzmann 
canonical ensemble can be computed as: 

, , /^di?^dP^aie-/5^(^"'^") 

(c?) = — (1) 

where (3 = l/ksT, ks is the Boltzmann's constant and T the temperature of the system. 
A^ is the number of particles in the system, R^ and P^ their coordinates and momenta 
respectively, and "H the Hamiltonian describing the system. Thus the aim of a CG method- 
ology is to build an effective Hamiltonian T-i once a set of {{u.i}) is known from either 
experiment or simulation. A generalization of this approach to several such variables will be 
described later in this letter. Under such assumptions, the Hamiltonian, l-i is rewritten as 
1-L = Href + J2f=i ^i^^i where the expectation values of K order parameters are known, Href 
is a reference Hamiltonian and A is a coupling constant. For simplicity in presentation, the 
theory is developed for a single order parameter, a = ai and K = 1. Then the expectation 
value can be written as: 



(^) = _ 1 _^ln /di?^dP^ ^~mr.AR-,P-Hx<r] 



The construction of a CG model can be thought of as equivalent to a minimization of the 
difference between the target (a) and the value [a^^) computed from the CG Hamiltonian 
1-L, namely: 

(a^«) - (a) = _A|iln|dil^dP^e-^[«^^/(^"'^")+^('^-<'^»l| 

= (3) 

where the condition A = was dropped. This equation should be regarded as a constitutive 
equation for CG models. Its meaning is straightforward: once Href is defined, the aim of 
the method is to find an optimal value for A that can best approximate the required (a). 
Another way to look at the problem is to recast Eq.Q in the following manner- If the 
expression inside the braces in Eq.Q is considered to be a free energy function, F, then Eq. 
(I3D has the form |f = 0. 



This form of Eq.([3]) suggests that the construction of a model representation can be 
pictured as follows. The observed order parameter average can be thought of as being 
induced by the perturbation of a reference system given by the Hamiltonian Tiref- This 
external perturbation is "conjugate" to the order parameter function, and finding a solution 
to Eq. ([2]) is equivalent to tuning the external perturbation of the reference system so as to 
obtain the target order parameter ensemble average. Thus, in order to obtain the required 
response a of the system, A is tuned. In the spirit of linear response theory, assuming that 
the perturbation is small, it is possible to write a Taylor expansion of Eq.(|3]) truncated at 
first order: 

+ OiXl) (4) 



{a}- {a}[0] = Ai 



9^ A=0 

where it is assumed that A = Ai such that (cr)[Ai] c^ (a). Eq.(jl]) can be solved iteratively. 
When the order parameter a is chosen to correspond to the pair distribution function, 
g{r) of a system, Eq.(j4]) becomes equivalent to the reverse monte carlo method. |j, l5| Thus, 
the reverse monte carlo approach is a natural consequence of the constitutive equation, 
Eq.([3]), when this equation is used to obtain a microscopic CG model that reproduces given 
structural expectation values as represented by the pair distribution functions. Further 
elaboration of this method can be found in the literature. [4, |5| 

Covariance matrices are a particular and important example of such correlation functions, 
and the y a re of critical importance for dimensionality reduction and data analysis in several 
fields Sl, [ll| . Thus, in general, if observations are made regarding the values of some function 
of the coordinates, {4>i{r^)}, then correlations (^[0j,0j]) between such functions can be 
written as: 

<^[0„0,] = 0,(r--)0,.(r--) (5) 

Thus, the goal becomes one wherein these expectation values have to be reproduced by the 
coarse grained potential. In the theory presented here, the CG potential energy for this 
system becomes: 

v, = T.Mnn<p^AJ] (6) 

Thus, Eq. ([3]) requires the construction of potentials Vtj (f") such that they satisfy the equa- 
tion: 






exp{-l3V^}dr^] =(^[0„0^.]) (7) 
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An approximate solution of this equation can be obtained through the iterative reverse 
monte carlo technique described previously (Eq. H]). An important special case involves the 
construction of CG hamiltonians that reproduce the covariance matrix of displacements, 
(f)i{r) = fi — (r*j), of atoms or sites labeled i, generated by a molecular ensemble. Such 
a covariance matrix is used to characterize conformational substates, and their collective 
fluctuation modes in proteins and other biomolecules. We now demonstrate that this ap- 
proach naturally leads to the construction of double and multiple well network models (9|, [l2| 
that have been proposed under phenomenological considerations, to describe mechanisms 
by which proteins undergo large scale conformational change. For simplicity, consider a 
situation wherein a protein structure is found to exist in two distinct conformations, which 
can be distinguished from each other by the value of some macroscopic order parameter, C,- 
For each of these values of .^ a different covariance matrix can be constructed. Then, one can 
construct a potential that reproduces the observed or simulated set of position covariances, 
^^[ri,fj] , or as a special case, root mean square fluctuations (B-factors) estimated from 
experiment, for each intermediate conformation (i.e each value of the order parameters, ^). 
Thus, the coarse-grained potential that reproduces these covariances becomes: 

Vr\Q = j:^^A^W,-Sr, (8) 

Experimental B-factors of proteins, that are proportional to mean square fluctuations of 
residues (i.e diagonal elements of the protein backbone covariance matrix) are often directly 
measured. Thus, a special case of such models is to construct CG Hamiltonians such that 
they reproduce experimentally measured B-factors of protein residues. An additional con- 
dition is imposed that the CG model of the protein structures retain the same topology as 



that observed in the original protein conformations 13J . With these conditions, if only two 



distinct conformational substates of the protein, characterized by a single order parameter 



^ G [0, 1] are considered, the "double well" network models [l2| used in the literature are 
recovered; Namely: 

v^A^]=VoT,,[0];^e [0,0.5] 

= ViT,,[l];^>0.5 (9) 

Here, vq and vi are constant and Tij[$^] is the connectivity matrix of the conformational 
substate C,, defined as Fjj = 1 if amino acid backbone atoms, i and j are within a prescribed 



cutoff, and zero otlierwise.|l3| Tlius, tliis tlieory provides a rigorous statistical mecfianical 
framework for tlie derivation of coarse-grained models that can be used to describe protein 
conformational transitions. 

It has been demonstrated that the theory presented can be used to construct microscopic 
statistical mechanical ensembles which reproduce low resolution information in the form of 
given expectation values of a few discrete order parameters. It will now be shown that the 
same theory can be utilized to perform the opposite and important process of reduction, 
wherein molecular coarse grained models are constructed from an equilibrium statistical 
mechanical ensemble of a complex system. This requires the application of the theory to 
the reproduction of continuous order parameter expectation values, since such reduction 
processes generally involve the construction of reduced many particle distribution functions, 
rather than discrete expectation values. For this purpose, the fundamental equation, Eq.([3]), 
is still valid once the following transformation is applied: 



\a-^ J d7^^v,(7^^)a(7^^|7^^) (lo) 

Correspondingly, all variations in the coupling constants. A, become functional derivatives: 

— ^ ^ (ID 

The "external potential" Vai'R-'^), Eq. ( jTOj) . will be used with the reference Hamiltonian, 
7iref{P^) = Z]/Li 2m"' where Mj is the mass associated with the /-th particle (site). In 
Eq.f lTOl) the symbol TZ^ is employed instead of i?^ in order to differentiate between the 
treatment of inverse methods and methods for reducing the degrees of freedom. In the 
latter the atomistic molecular ensemble of the system is known, and the target of this 
methodology is to build a new coarse grained model associated with a simplified topology 
of the original system, which is consistent with the original molecular ensemble |6|. To 
construct such a model, a mapping of the position of the atomistic coordinates r" to CG 
sites in the simplified topology, TZ^ , is defined through a mapping operator M^m{t^)- This 
mapping operator relates the two sets of coordinates through 

n^ = M^{r^). (12) 

An example of such a mapping is one where the coordinates r of the atoms are mapped 
onto their center of mass. The corresponding CG potential defined by the right side of 

6 



Eq. (fTOj) has a simple physical interpretation. The integral computes the potential energy 
of a particular configuration (or "snapshot") of the system. Each such configuration is 
represented by the variable TZ^ . Since the order parameter, (j{TZ^\R.'^) is a continuous 
function defined over the coarse grained space, the potential energy due to the atoms at 
positions TZ^ is integrated over all the spatial degrees of freedom in Eq. flTUl) . The energy 
of each such snapshot can then be written as the sum over interactions between each atom 
at the position TZ ^ and a corresponding external field Va{TZ ^). The external field couples 
to the atoms through the order parameter a{TZ^\R.^): for example, an electric field will 
interact with the charge of atoms, whereas a magnetic field will affect the corresponding 
spin state of each atom. This potential energy associated with a particular configuration (or 
"snapshot") TZ^ is the "mean-field". 

As with models constructed using inverse methods, the problem of building a CG model 
consists in solving Eq. (|3]). If Z is the partition function associated with the CG Hamiltonian 

•H(i?^,P^) = 'HrefiP'') + /d^'^V.(^'^)a(7^'^|7^^) 

= nref{P'') + npot{n^) (is) 

then Eq.(l3]) can be recast as: 

= {Hn^))AA. (14) 

Here Z^°"'^[a] = JdTZ^e'^^^^ \ This equation can be solved variationally. A specific 
choice for (o"(7^o^)) that has been used in literature is to minimize the difference between 
the many particle reduced distribution function, (cr(7?.^))^A = {S[7Zq — M^(f")]), and the 
many particle distribution function in the simplified topology, {o'{1Zq))cg = {^[^o ~ ^ ^]) 
for V7^^ G TZ^ . For this specific selection the effective potential can be written as: 

. yconf _ 



V7^o(^o) = -^B^ In 



^ _ A/f^(^\-\.-fi^{^) 



dr="5[7^;; - M^if-y^e 



+ ln 



yConf 



(15) 



where 2;™"/ = ^ df~^^-P'^{^) and u(P') is the atomistic potential energy. Observing that the 
potential is defined upto a constant, and that /c^T" In ^1°"^ = const, this equation corresponds 
to the fundamental result underlying the MS-CG method 6|. Elaborations of this method 
can be found in the cited literature. 



Both the inverse methods and the CG formulations derived here as specific implemen- 
tations of Eq. ([3]) involve in practice, a minimization of the difference in the expectation 
values of an order parameter when evaluated over a CG ensemble relative to an ensemble 
average measured either from experiment or simulation. A natural means by which such a 
minimization can be achieved is to perform a linear least squares fit of the unknown (coarse- 
grained) average to the target value. Such a fitting procedure offers a simple but powerful 
numerical implementation of the constitutive equation, Eq. ([3]) and has been extensively used 
in the context of constructing molecular CG models. 6| 

Eq. flT^ can also be interpreted in terms of the loss of information relative to the original 
atomistic ensemble, upon performing the reduction procedure. This equation states that 
the Boltzmann factor associated with a specific CG conformation TZ^ should be equal to 
the number of microscopic conformations that can map to it, divided the total number 
of microscopic states. This interpretation of Eq. (fT5|) suggests an alternative approach to 
construct CG representations through the minimization of the relative entropy between the 
atomistic and lower resolution CG system. Such a formalism has been recently employed 
|l4| to derive CG potentials. Within the relative entropy formalism it can be shown that 
the minimal entropic contribution due to this loss of information, Smap, is given by: 

Smap = -kB\nl^^Jdf^6[n^ - M^(r--)]e-^"(^^)| (16) 

Thus, the entropic contribution is (upto a constant) nothing but the logarithm of the CG 
potential energy in Eq.( TT5|) . It should be noted that relative entropy arguments can be 
used to recover the Gibbs-Bogoliubov-Feynman thermodynamic inequality that is of crucial 



•3 



importance in equilibrium statistical mechanics. 

As a final example of the generality of the approach presented here, it is easy to see 
that if the ensemble averaged single particle density of a system, p(r) = -^ J27=i {^{^~ ^i)) 
is chosen to be the desired target order parameter average, then Eq.([3]) implies that this 
requires the construction of a single particle potential v{r) = v[p]{f) that must generate 
a Boltzmann ensemble which reproduces this average density. The potential derived as a 
result of this requirement is thus a functional of the single particle density, and the model 
formulation that results corresponds to the well known classical density functional theory of 



classical statistical mechanics. 
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In conclusion, it has been shown that many commonly used methods to analyse classical 



statistical ensembles are specific realizations of a generalized mean field theory that can be 
used to model these ensembles. The generalized formalism proposed here also suggests that 
equilibrium coarse graining methods that are routinely employed in the study of complex 
molecular systems are closely linked to inverse problems involving these systems, where low 
resolution information from experiment or simulation is used to reconstruct a microscopic 
description of the underlying physical system. It has also been shown that both atomistic 
density functional theories and lower resolution mean field coarse grained models can be 
derived from this formalism through an appropriate choice of the target properties of in- 
terest. Thus, the generalized mean field theory discussed in this work offers a simple and 
general framework for both the construction of sophisticated coarse grained approaches to 
modeling the statistical mechanical properties of complex molecular systems, as well as the 
reconstruction of microscopic models based on macroscopic experimental observations of 
their ensemble averaged properties. 

We thank Hans Andersen and Gregory Voth for valuable comments and discussions. 
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